Vortex Line Ordering in the Driven Three-Dimensional Vortex Glass 
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Resistively-shunted-junction dynamics is applied to the three dimensional uniformly frustrated 
XY model with randomly perturbed couplings, as a model for driven steady states in a type-II 
superconductor with quenched point pinning. For a disorder strength p strong enough to produce a 
vortex glass in equilibrium, we map the phase diagram as a function of temperature T and uniform 
driving current I. Using finite size analysis and averaging over many realizations of quenched 
randomness we find a first-order melting T m (J) from a vortex line smectic to an anisotropic liquid. 
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Ordering and phase transitions in driven steady states 
far from equilibrium remains a topic of considerable 
general interest. In particular, the spatial ordering of 
driven vortex lines in a type-II superconductor with ran- 
dom point pinning has received considerable theoreti- 
cal LU, 

H, 0, 0] and experimental [H attention. Origi- 
nally, Koshelev and Vinokur proposed [l[ that a mov- 
ing steady state would average over quenched random- 
ness, and that a system which was disordered in equi- 
librium could reform into an ordered vortex lattice when 
driven. Later, Giamarchi and Le Doussal [2J argued that 
this state would be a "moving Bragg glass," with alge- 
braically decaying translational correlations both parallel 
and transverse to the driving force. Balents, Marchetti 
and Radzihovsky Q then argued that the moving Bragg 
glass would be unstable to dislocations which decouple 
the planes of vortex lines moving parallel to the drive, 
resulting in a smectic ordering. Other theoretical works 
[H have supported one or more of these scenarios. 

While numerous simulations have studied this prob- 
lem for point vortices in two dimensions [||, few works 
have treated three dimensional (3D) vortex lines at finite 
temperature @, 0, 0, E3] ■ We present here the first 3D 
simulations to include both a systematic study of finite 
size effects, as well as averaging over many realizations of 
the quenched randomness Such considerations are 

necessary for an unambiguous determination of ordering 
in driven steady states. We map the phase diagram as a 
function of temperature T and driving current / for the 
strongly random case, and study the nature of ordering 
just below and above the melting transition T m (I). 

Following Refs. [1,(1, [ijj we use the 3D XY modeljll 
with resistively shunted junction (RSJ) dynamics [13j to 
model our system. The Hamiltonian is given by, 



H[0(Ti)] = -J2 J ^ cos ( e ( r *) - + A) - A ifl ) (1) 



function on site Ti of a cubic L x x L y x L z grid of sites 
with bonds in directions \i = x, y, z. The circulation of 
Ai^ around any plaquette of the grid is fixed and equal to 
27r/ with / the density of applied magnetic flux quanta 
through that plaquette. We use a uniform value / = 1/12 
oriented in the z direction. The resulting field induced 
vortex line density, also equal to /, forms a vortex lattice 
in the equilibrium ground state of the pure (disorder-free) 
system, with basis vectors ai = Ay and a.2 = ix + 2y. To 
model quenched point randomness we use couplings 



distributed, 



"P £ i/i)> where are uncorrelated, uniformly 



random variables with (e» 



0, 



1. 



We use Jfj_ = J 



in the xy plane, and J z — J±/A0 to en- 
hance vortex line fluctuations along the z direction. The 
disorder strength is controlled by the parameter p. In our 
earlier work 
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15[ we showed that above a critical p c 
the vortex line lattice becomes unstable to a vortex glass 
at low temperatures in equilibrium. Here we consider the 
strongly random limit p — 0.15 > p c « 0.14. 

We use RSJ dynamics, with equation of motion as in 
9fl, however with fluctuating twist boundary conditions 
13j | in all directions. We apply a uniform current I = Ix, 
resulting in a force ixlon the vortex lines driving them 
in the y direction. We use a second order Runge-Kutta 
integration method with dimensionless time step [9| At = 
0.1, and typically 2.6 x 10 6 time steps per simulation 
run, resulting in a net displacement of the vortex line 
center of mass of roughly 48, 000 grid spacings in the 
ordered phase. Depending on system size, up to 3/4 of 
these steps may be discarded for equilibration. When 
probing behavior at a specific point in the T — I plane, we 
usually take the pure system ground state as our initial 
configuration. For numerous test cases, however, we have 
started with a random initial state at high T and slowly 
cooled to the desired point. Except for a narrow region 
of hysteresis at the melting transition, we always find the 
same long time steady state for both initial conditions. 



where 0{ri) is the phase of the superconducting wave- 



To determine the structural order of our system we use 
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the vortex structure function, 



L-^e lh {, 2 (R + r)n 2 (R)> (2) 



fL x L y L z 



R,r 



where n z (r) is the vortex line density in the z direction 
at position r. We use (. . .) to denote averages over simu- 
lation time and [. . .] to denote averages over independent 
realizations of quenched randomness (typically 40 real- 
izations are used). We also consider the correlations, 



C(x,y,z) = _L-]rs(k) e 

L/ x L/ y L/ z 



(3) 



C ( , . /,,, . -. ) = _L- S(k)e-« k * x + k ^ (4) 



Using the appearance of sharp peaks in S^kx, k z = 0) 
to signal an ordered phase, in Fig. [T^, we present the T—I 
phase diagram for a 24 x 24 x 16 sized system for both 
the pure (p = 0) and random (p = 0.15) cases. We mea- 
sure T in units of J±, and / in units of Iq = 2eJ±/h. 
Crossing the phase boundary at any value of the current 
we find a discrete jump in energy, suggesting a first or- 
der melting transition T m {I). For the random case our 
phase boundary is from a single random realization only 
and is determined by two separate methods: by cooling 
in T at fixed / from the disordered phase (solid line); by 
increasing and decreasing I at fixed T from within the or- 
dered phase (dashed line). Assuming some hysteresis, as 
in a first order transition, the two methods give reason- 
able agreement. At fixed T, the random system orders 
only for a finite current interval I c i{T) < I < I c2 (T). 
For the pure case, our Tm{I) looks qualitatively simi- 
lar to that found in Ref.[9j for the weakly random case 
p < p c . The melting of the ordered phase upon increas- 
ing I > I C 2, which exists for both random and pure sys- 
tems, is contrary to theoretical predictions [H, 0, 0, |j] . 
Ref. suggested this to be a consequence of thermally 
excited vortex loops. We in fact find that the density 
of thermally excited loops is negligible everywhere in the 
ordered phase, but shows a finite jump to a noticeable 
amount as T m (I) is crossed; this jump is particularly 
large as / crosses I C 2(T). 

Consider first the disordered phase. Fig. []J> shows an 
intensity plot of lnS^kj^A^ = 0) for a single realiza- 
tion of a random 48 3 system at / = 0.48, T = 0.13, 
just above the maximum in T m {I). Two relatively faint 
peaks, S(K )/S(0) ~ 0.008, lie off the origin along the k x 
axis. Previous works 0, [Toj ] have interpreted such peaks 
as evidence for a smectic phase. However, as shown in 
Fig. [2}d, we find no significant increase of peak height 
5(Kq) with system size, thus indicating that this is an 
anisotropic liquid with only short ranged translational 
correlations. 

Next consider the ordered phase. Fig. [5^ shows an in- 
tensity plot of In S(k± , k z = 0) for a single realization of a 
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FIG. 1: (a) Phase diagram as function of temperature T and 
uniform driving current Ix (vortex lines move in y direction). 
Solid symbols are for a pure system, p = 0. Open symbols are 
for a system with quenched point pinning, p — 0.15; o is phase 
boundary obtained by cooling T from the disordered phase at 
fixed /, while □ is phase boundary obtained by increasing 
or decreasing / from the ordered phase at fixed T. Results 
are from a 24 x 24 x 16 size system, (b) Intensity plot of 
\n S(k±, k z = 0), for one particular random realization (p — 
0.15), in disordered phase at I = 0.48, T = 0.13 for a 48 x 
48 x 48 size system; k = is at the center of the figure. 



random 48 3 system at I — 0.48, T = 0.09, just below the 
maximum in T m (I). Now sharp peaks, iS , peak/S'(0) ~ 1, 
lie along the k x axis. In the pure system, the vortex 
lattice state has Bragg peaks at Kiq = (2n/3)x and 
K u = (2tt)(.t/6 + y/4), as labeled in Fig. Hi. In the 
random system a peak remains exactly at Kio; a sec- 
ond peak is always found with K y — K\\^ v but with 
K x — K\i lX — 0,±2ir/L y depending on the particular 
random realization. We therefore generalize our notation 
to denote by Ku the exact location of this second peak, 
with the understanding that the value of K\\ x may shift 
slightly between different random realizations. In Fig. [2)3 
we plot the disorder averaged height of these peaks, 
[S(Kio)] and [S(Kn)], vs. system volume V — L x L y L z , 
for several different system sizes, for both the pure and 
random cases. In both cases [<S(Kio)] scales linearly with 
V, indicating a sharp Bragg peak. This Bragg peak in- 
dicates that vortex lines are organized into specific yz 
planes with a periodic spacing of 3i between planes. Or- 
dering within and between planes is reflected in the scal- 
ing of [S(Kn)]. For the pure case, S(Ku) ~ V indicat- 
ing the long range translational order of a moving vor- 
tex line lattice. For the random case, however, [S(Kn)] 
grows less rapidly than V. The dashed line in Fig. [5b rep- 
resents a power law divergence of [S(Kn)] ~ p 2 / 3 ~ L 2 , 
however if we discard the smallest point at L = 12, the 
exponent of a power law fit in L becomes 1.65. It is thus 
unclear from Fig. [2Jd whether our biggest size system has 
reached the asymptotically large L limit. The sublinear 
growth of [5(Kn)] with V indicates that the system has 
less than long range translational order. However it does 
not indicate which direction(s) are less ordered than oth- 
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FIG. 2: (a) Intensity plot of In i5(kj_, k z — 0), for one particu- 
lar random realization, in ordered phase at / = 0.48, T = 0.09 
for a L x = L y = L z = 48 size system; k = is at the center 
of the figure, (b) Scaling of disorder averaged peaks heights 
vs. system volume V = L x L y L z . For the ordered state of (a) 
we plot [5(Kio)] and [S(Ku)j for pure (p = 0) systems of 
size 24 x 24 x L (squares) and L x L x 24 (triangles), and for 
systems with quenched point pinning (p — 0.15) of size L 3 , 
L = 12,24,36,48,60 (circles). We also plot [S(K )] for the 
p — 0.15 disordered state at T = 0.13 of Fig. fjj) for systems 
of size L 3 with L = 24, 36, 48, 60 (crosses). 

ers. For this we will consider the shape of the peak at 
Kn and the real space correlations of Eqs. and (T5|). 

The sharpness of the peak at Kn in the k y direction 
suggests that vortex lines are periodically ordered along 
their direction of motion y within each yz plane. That 
this peak appears broad in the k x direction suggests that 
the yz planes have only short range correlations between 
them. We thus refer to the yz planes containing the vor- 
tex lines as smectic planes. We now consider the real 
space correlations of Eqs. ([3]) and Fig. shows an 
intensity plot of C(x,y, 2 = 0) corresponding to S(k) of 
Fig. [5^. We clearly see the structure proposed above: 
vortices lie in periodically spaced yz planes with separa- 
tion 3x; within a given plane, i.e. x — 0, vortices are pe- 
riodic with separation Ay] at finite transverse separation, 
i.e. |x| > 0, the variation of C(x,y, 0) with y decreases, 
until at large \x\ it is almost uniform in y, thus indicating 
short ranged correlations between the smectic planes. To 
quantify this, we consider the correlation Eq. Q evalu- 
ated at k y = Kn t y giving the periodicity within a given 
smectic plane. Since the slightly different values of Ku x 
for the different random realizations give rise to differ- 
ent complex phase shifts in C{x,Kn^ y ,z = 0), we dis- 
order average the absolute value. In Fig. [3J3 we plot 
[\C(x, Ku y y,0)\}/Co vs. x, for different system sizes L 3 . 
We normalize by Co = [\C(0, Ku iV , 0)|] to better com- 
pare different system sizes. We see a clear exponential 
decay to zero; solid lines are fits to a periodic exponen- 
tial giving correlation lengths in the range £ x ~ 5 — 6.5. 
In recent works [tj [l(| pictures similar to Fig. [2^ were 
identified as a "moving Bragg glass." However the short 
range correlations of Fig. [3Jd clearly show our ordered 
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FIG. 3: (a) Intensity plot of real space correlation, C(x,y, z — 
0), for one particular random realization, in ordered phase at 
/ = 0.48, T = 0.09, p = 0.15, for a L x = L y = L z = 48 size 
system; r = is at the center of the figure, (b) Corresponding 
disorder averaged correlation, [\C(x, Kn, y , z — 0)|] vs. x for 
systems of size L x L x L. 



phase to be a smectic rather than a Bragg glass. 

Next we consider correlations parallel to the direction 
of motion y, within a given smectic plane. Fig. 2£i shows 
[C(0, y, 0)] vs. y, for various system sizes, for the same 
I = 0.48, T = 0.09 as in Figs. and H Since we found 
l; x ~ 6 is fairly small, we include also systems of size 
36 x L x L in order to obtain larger values of L y . Fig. Q] 
suggests decay to a periodic oscillation, with only a small 
finite size effect when y ~ L y /2. The large y limit can be 
described by the three envelop functions, [C max (Lj / )] = 
min m [C(0,4m,0)], [C min (L y )} = max m [C(0, 4m + 2, 0)], 
and [C m id(L y )] = max m [C(0,4m ± 1,0)], for m integer, 
which we plot vs. L y in Fig. UJd. For long range transla- 
tional order, these should converge to different constants 
as L y — > oo. For algebraic or short range order, these 
should all converge to the value 1 /4 characterizing a uni- 
form vortex line density. We find that both exponential 
decay to long range order (red solid line), and an alge- 
braic decay to 1/4 (blue dashed line), give equally good 
fits. For long range order, the fit to [C max ] gives a decay 
length ^ y ~ 25, while the algebraic fit to [C max ] gives a 
power law exponent of ~ 0.2. 

Finally we consider the correlations along z, parallel to 
the applied magnetic field, within a given smectic plane. 
Fig. [SJi shows [(7(0, 0, z)] — 1/4 vs. z, for various system 
sizes, for the same I = 0.48, T = 0.09 as before. The 
subtracted value 1/4 is the uniform vortex density in the 
smectic plane, in the absence of any ordering. For the 
largest system size L y = 96 we find an exponential decay 
to zero with a decay length £ z ~ 9. Smaller sizes show 
a decay to a small finite constant that decreases with 
increasing L yi consistent with our earlier observation in 
Fig.[2jD that [S(Ku)] has not yet reached the asymptotic 
large L limit for L < 60. In Fig. [SJd we show an in- 
stantaneous configuration of vortex lines in a particular 
smectic plane. At fixed z, the lines appear periodically 
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FIG. 4: (a) Disorder averaged real space correlation, [C(x = 
0,y,z — 0)] vs. y for various system sizes, in the ordered 
phase at I = 0.48, T = 0.09, p = 0.15. (b) Corresponding 
limiting values of maximum, middle, and minimum envelop of 
[C(x = 0,y,z = 0)] (see text) vs. L y for various system sizes; 
red solid line is a fit to an exponential decay to a constant, 
blue dashed line is a fit to an algebraic decay to 1 /4. 
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FIG. 5: (a) Disorder averaged real space correlation, [C(x = 
0, y = 0, z)] — 1/4 vs. z for various system sizes, in the ordered 
phase at I = 0.48, T = 0.09, p = 0.15. (b) Instantaneous 
vortex line configuration in a particular smectic plane of a 
particular random realization, in the ordered phase at the 
same parameters as (a), for a system of size 36 x 96 x 96. 



spaced along y, in agreement with Fig. |4j However, trac- 
ing any particular line along z, we see that it wanders 
an amount comparable to the inter-line spacing, consis- 
tent with the decay seen in Fig. In some planes, such 
as the one shown in Fig. [SJd, lines can have a net tilt, 
by closing onto one of their neighbors under the peri- 
odic boundary condition along z. Each such tilted plane 
must be compensated by another plane that tilts in the 
opposite direction, so that the net vorticity in the y di- 
rection vanishes. Looking at individual configurations, 
we find a strong relation between the correlations along 
z within a given smectic plane, and the correlations along 
x between different smectic planes. When the lines in a 
smectic plane are tilted, as in Fig. [5Jd, or otherwise have 
a transverse wandering comparable to the average spac- 
ing between lines, we find that this plane decouples from 



its neighbors, moving either diffusively with respect to 
its neighbors, or even with a slightly different average 
speed. Planes whose lines have small transverse wander- 
ings, remain strongly correlated with their neighbors. It 
is necessary to have only a few such decoupled smectic 
planes in order to have short range average correlations 
along x, as in Fig. [3Jd. 

Our analysis of the ordered phase has been for a point 
just below the peak in the first order melting curve T m (I). 
Establishing the nature of correlations throughout the 
ordered phase would require similar finite-size, random- 
realization-averaged, analyses elsewhere in this region. 
Our preliminary investigations suggest that as either T or 
/ is decreased, finite correlation lengths grow and become 
too large to make such studies feasible at present. We 
thus cannot rule out the possibility of a more ordered 
"moving Bragg glass" phase at lower temperatures. 
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